#include "collision.h"
#include "primitive.h"

using namespace Collision;

#include <glm/glm.hpp>
#include <glm/gtx/norm.hpp>

#include <QDebug>


glm::vec4 closestPointOBB(Box &cube, glm::vec4 dest, glm::mat4 &objToLocal)
{
    glm::vec4 temp_points[8];
    for(int i = 0; i < 8; i++) {
        temp_points[i] = objToLocal * cube.m_points[i];
    }
    glm::vec4 temp_origin = objToLocal*cube.m_origin ;
    glm::vec4 distance = dest - temp_origin;

    glm::vec4 result = temp_origin;

    glm::vec4 box_Axis[3];
    float half_side[3];

    box_Axis[0] = glm::normalize(temp_points[1] - temp_points[0]);
    half_side[0] = std::sqrt(glm::gtx::norm::length2(temp_points[1] - temp_points[0]))*0.5;
    box_Axis[1] = glm::normalize(temp_points[3] - temp_points[0]);
    half_side[1] = std::sqrt(glm::gtx::norm::length2(temp_points[3] - temp_points[0]))*0.5;
    box_Axis[2] = glm::normalize(temp_points[4] - temp_points[0]);
    half_side[2] = std::sqrt(glm::gtx::norm::length2(temp_points[4] - temp_points[0]))*0.5;

    //qDebug() << "Distance: " << distance.x << " " << distance.y << " " << distance.z;

    // For each box axis
    for ( int i = 0; i < 3; i++ ) {
        float axis_projection = glm::dot(distance, box_Axis[i]);
        //qDebug() << "Box Axis " << i << ": " << box_Axis[i].x << " "  << box_Axis[i].y << " " << box_Axis[i].z << " ";
        //qDebug() << "Axis " << i <<  " projection: " << axis_projection;
        //qDebug() << "Half-side " << i << ": " << half_side[i];

        // Clamp to the box extents
        if ( axis_projection > half_side[i]) {
            axis_projection = half_side[i];
        }
        if ( axis_projection < -half_side[i]) {
            axis_projection = -half_side[i];
        }
        result += axis_projection*box_Axis[i];
    }
    return result;
    /*
    glm::vec4 temp_points[8];
    for(int i = 0; i < 8; i++) {
        temp_points[i] = objToLocal * temp_points[i];
    }
    glm::vec4 temp_origin = objToLocal*cube.m_origin ;
    glm::vec4 distance = dest - temp_origin;

    glm::vec4 result = temp_origin;

    glm::vec4 box_Axis[3];
    float half_side[3];

    box_Axis[0] = glm::normalize(temp_points[1] - temp_points[0]);
    half_side[0] = std::sqrt(glm::gtx::norm::length2(temp_points[1] - temp_points[0]))*0.5;
    box_Axis[1] = glm::normalize(temp_points[3] - temp_points[0]);
    half_side[1] = std::sqrt(glm::gtx::norm::length2(temp_points[3] - temp_points[0]))*0.5;
    box_Axis[2] = glm::normalize(temp_points[4] - temp_points[0]);
    half_side[2] = std::sqrt(glm::gtx::norm::length2(temp_points[4] - temp_points[0]))*0.5;

    //qDebug() << "Distance: " << distance.x << " " << distance.y << " " << distance.z;

    // For each box axis
    for ( int i = 0; i < 3; i++ ) {
        float axis_projection = glm::dot(distance, box_Axis[i]);
        //qDebug() << "Box Axis " << i << ": " << box_Axis[i].x << " "  << box_Axis[i].y << " " << box_Axis[i].z << " ";
        //qDebug() << "Axis " << i <<  " projection: " << axis_projection;
        //qDebug() << "Half-side " << i << ": " << half_side[i];

        // Clamp to the box extents
        if ( axis_projection > half_side[i]) {
            axis_projection = half_side[i];
        }
        if ( axis_projection < -half_side[i]) {
            axis_projection = -half_side[i];
        }
        result += axis_projection*box_Axis[i];
    }
    return result;
    */
}

CollisionState Collision::sphere_cube_collision(Sphere &sphere, Box &cube)
{
    glm::mat4 cube_to_sphere = sphere.m_frame->world_to_local()*cube.m_frame->local_to_world();
    //glm::mat4 cube_to_sphere = sphere.getToLocalMatrix()*cube.getToWorldMatrix();

//    qDebug() << "C0: " << cube_to_sphere[0].x << " "  << cube_to_sphere[0].y << " "  << cube_to_sphere[0].z << " "  << cube_to_sphere[0].w;
//    qDebug() << "C1: " << cube_to_sphere[1].x << " "  << cube_to_sphere[1].y << " "  << cube_to_sphere[1].z << " "  << cube_to_sphere[1].w;
//    qDebug() << "C3: " << cube_to_sphere[2].x << " "  << cube_to_sphere[2].y << " "  << cube_to_sphere[2].z << " "  << cube_to_sphere[2].w;
//    qDebug() << "C4: " << cube_to_sphere[3].x << " "  << cube_to_sphere[3].y << " "  << cube_to_sphere[3].z << " "  << cube_to_sphere[3].w;

    glm::vec4 closest = closestPointOBB( cube, sphere.m_origin, cube_to_sphere );
    //qDebug() << "Box Data: " << cube.m_origin.x << ", " << cube.m_origin.y << ", " << cube.m_origin.z << " and Side: " << cube.m_side;
    //qDebug() << "Closest Point: " << closest.x << ", " << closest.y << ", " << closest.z;

    if ( (glm::gtx::norm::distance2( closest, sphere.m_origin ) > sphere.m_radius*sphere.m_radius)  ) {
        //qDebug() << "Distance " << glm::gtx::norm::distance2( closest, sphere.m_origin ) << " - radius: " << sphere.m_radius*sphere.m_radius ;
        return Outside;
    }

    //int inside_pts = 0;
    for (int i = 0; i < 8; i++ ) {

        if (!sphere.inside(cube_to_sphere*cube.m_points[i])){///*sphere.inside(sphere.getToLocalMatrix()*cube.getToWorldMatrix()*cube.m_points[i])*/) {
            return Partially;
        }
    }
    return Inside;
    /*
    //qDebug() << "Number of points inside the sphere: " << inside_pts;
    if ( inside_pts == 8 ) {
        return Inside;  // Box Inside sphere
    } else {
        return Partially;
    }*/
}


// box - the obb that we are testing agains the - its coordinate system will be used as reference
// box_node - the tree box
CollisionState Collision::obb_obb_collision( Box & box, Box &box_node )
{
    CollisionState ret_val = Outside;
    glm::mat4 node_to_box = box.m_frame->world_to_local()*box_node.m_frame->local_to_world();
    //glm::mat4 cube_to_sphere = sphere.m_frame->world_to_local()*cube.m_frame->local_to_world();
    qDebug() << "C0: " << node_to_box[0].x << " "  << node_to_box[0].y << " "  << node_to_box[0].z << " "  << node_to_box[0].w;
    qDebug() << "C1: " << node_to_box[1].x << " "  << node_to_box[1].y << " "  << node_to_box[1].z << " "  << node_to_box[1].w;
    qDebug() << "C3: " << node_to_box[2].x << " "  << node_to_box[2].y << " "  << node_to_box[2].z << " "  << node_to_box[2].w;
    qDebug() << "C4: " << node_to_box[3].x << " "  << node_to_box[3].y << " "  << node_to_box[3].z << " "  << node_to_box[3].w;


    glm::vec4 temp_points[8];
    for(int i = 0; i < 8; i++) {
        temp_points[i] = node_to_box * box_node.m_points[i];
    }
    glm::vec4 temp_origin = node_to_box*box_node.m_origin ;


    // 0, 1, 2 are lines from GToA matrix
    glm::vec4 GToA[3];
    GToA[0] = glm::normalize(box.m_points[1] - box.m_points[0]);
    GToA[1] = glm::normalize(box.m_points[3] - box.m_points[0]);
    GToA[2] = glm::normalize(box.m_points[4] - box.m_points[0]);

    glm::vec4 extentsA(0.5f*std::sqrt(glm::gtx::norm::distance2(box.m_points[1], box.m_points[0])),
                       0.5f*std::sqrt(glm::gtx::norm::distance2(box.m_points[3], box.m_points[0])),
                       0.5f*std::sqrt(glm::gtx::norm::distance2(box.m_points[4], box.m_points[0])),
                       0.0);
    glm::vec4 GToB[] = {glm::normalize(temp_points[1] - temp_points[0]),
                         glm::normalize(temp_points[3] - temp_points[0]),
                         glm::normalize(temp_points[4] - temp_points[0])};
    glm::vec4 extentsB(0.5f*std::sqrt(glm::gtx::norm::distance2(temp_points[1], temp_points[0])),
                       0.5f*std::sqrt(glm::gtx::norm::distance2(temp_points[3], temp_points[0])),
                       0.5f*std::sqrt(glm::gtx::norm::distance2(temp_points[4], temp_points[0])),
                       0.0);

    glm::vec4 translation = temp_origin - box.m_origin;

    glm::vec4 translationInA( glm::dot(translation,GToA[0]),
                              glm::dot(translation,GToA[1]),
                              glm::dot(translation,GToA[2]),
                              0.0f );

    glm::mat3x3 R;  // B to A

    for(int i = 0; i < 3; i++) {
        for(int k = 0; k < 3; k++)
        {
            R[k][i] = glm::dot(GToA[i], GToB[k]);
        }
    }
    float ra, rb, t;
    long i, k;
    // Separatting Axis
    // A basis vector
    for( int i = 0; i < 3; i++ ) {
        ra = extentsA[i];

        rb = extentsB[0]*std::abs(R[0][i]) + extentsB[1]*std::abs(R[1][i]) + extentsB[2]*std::abs(R[2][i]);

        t = std::abs( translationInA[i] );
        if ( t > ra + rb ) {
            return Outside;
        }
    }
    // B basis vector
    for( int k = 0; k < 3; k++ ) {
        ra = extentsA[0]*std::abs(R[k][0]) + extentsA[1]*std::abs(R[k][1]) + extentsA[2]*std::abs(R[k][2]);
        rb = extentsB[k];

        t = std::abs( translationInA[0]*R[k][0] + translationInA[1]*R[k][1] + translationInA[2]*R[k][2] );

        if( t > ra + rb )
            return Outside;
    }

    // 9 cross products - After testing agains each obb axis

    // L = A0 x B0
    ra =
    extentsA[1]*std::abs(R[0][2]) + extentsA[2]*std::abs(R[0][1]);

    rb =
    extentsB[1]*std::abs(R[2][0]) + extentsB[2]*std::abs(R[1][0]);

    t =
    std::abs( translationInA[2]*R[0][1] -
         translationInA[1]*R[0][2] );

    if( t > ra + rb )
        return Outside;

    //L = A0 x B1
    ra =
    extentsA[1]*std::abs(R[1][2]) + extentsA[2]*std::abs(R[1][1]);

    rb =
    extentsB[0]*std::abs(R[2][0]) + extentsB[2]*std::abs(R[0][0]);

    t =
    std::abs( translationInA[2]*R[1][1] -
         translationInA[1]*R[1][2] );

    if( t > ra + rb )
        return Outside;

    //L = A0 x B2
    ra =
    extentsA[1]*std::abs(R[2][2]) + extentsA[2]*std::abs(R[2][1]);

    rb =
    extentsB[0]*std::abs(R[1][0]) + extentsB[1]*std::abs(R[0][0]);

    t =
    std::abs( translationInA[2]*R[2][1] -
         translationInA[1]*R[2][2] );

    if( t > ra + rb )
        return Outside;

    //L = A1 x B0
    ra =
    extentsA[0]*std::abs(R[0][2]) + extentsA[2]*std::abs(R[0][0]);

    rb =
    extentsB[1]*std::abs(R[2][1]) + extentsB[2]*std::abs(R[1][1]);

    t =
    std::abs( translationInA[0]*R[0][2] -
         translationInA[2]*R[0][0] );

    if( t > ra + rb )
        return Outside;

    //L = A1 x B1
    ra =
    extentsA[0]*std::abs(R[1][2]) + extentsA[2]*std::abs(R[1][0]);

    rb =
    extentsB[0]*std::abs(R[2][1]) + extentsB[2]*std::abs(R[0][1]);

    t =
    std::abs( translationInA[0]*R[1][2] -
         translationInA[2]*R[1][0] );

    if( t > ra + rb )
        return Outside;

    //L = A1 x B2
    ra =
    extentsA[0]*std::abs(R[2][2]) + extentsA[2]*std::abs(R[2][0]);

    rb =
    extentsB[0]*std::abs(R[1][1]) + extentsB[1]*std::abs(R[0][1]);

    t =
    std::abs( translationInA[0]*R[2][2] -
         translationInA[2]*R[2][0] );

    if( t > ra + rb )
        return Outside;

    //L = A2 x B0
    ra =
    extentsA[0]*std::abs(R[0][1]) + extentsA[1]*std::abs(R[0][0]);

    rb =
    extentsB[1]*std::abs(R[2][2]) + extentsB[2]*std::abs(R[1][2]);

    t =
    std::abs( translationInA[1]*R[0][0] -
         translationInA[0]*R[0][1] );

    if( t > ra + rb )
        return Outside;

    //L = A2 x B1
    ra =
    extentsA[0]*std::abs(R[1][1]) + extentsA[1]*std::abs(R[1][0]);

    rb =
    extentsB[0] *std::abs(R[2][2]) + extentsB[2]*std::abs(R[0][2]);

    t =
    std::abs( translationInA[1]*R[1][0] -
         translationInA[0]*R[1][1] );

    if( t > ra + rb )
        return Outside;

    //L = A2 x B2
    ra =
    extentsA[0]*std::abs(R[2][1]) + extentsA[1]*std::abs(R[2][0]);

    rb =
    extentsB[0]*std::abs(R[1][2]) + extentsB[1]*std::abs(R[0][2]);

    t =
    std::abs( translationInA[1]*R[2][0] -
         translationInA[0]*R[2][1] );

    if( t > ra + rb )
        return Outside;

    /*no separating axis found,
     the two boxes overlap */

    for (int i = 0; i < 8; i++){
        if(!box.inside(temp_points[i]))
            return Partially;
    }
    return Inside;



//    CollisionState ret_val = Outside;

//    glm::mat4 cube_to_sphere = sphere.m_frame->world_to_local()*cube.m_frame->local_to_world();

//    // 0, 1, 2 are lines from GToA matrix
//    glm::vec4 GToA[3];
//    GToA[0] = glm::normalize(box.m_points[1] - box.m_points[0]);
//    GToA[1] = glm::normalize(box.m_points[3] - box.m_points[0]);
//    GToA[2] = glm::normalize(box.m_points[4] - box.m_points[0]);

//    glm::vec4 extentsA(0.5f*std::sqrt(glm::gtx::norm::distance2(box.m_points[1], box.m_points[0])),
//                       0.5f*std::sqrt(glm::gtx::norm::distance2(box.m_points[3], box.m_points[0])),
//                       0.5f*std::sqrt(glm::gtx::norm::distance2(box.m_points[4], box.m_points[0])),
//                       0.0);
//    glm::vec4 GToB[] = {glm::normalize(box_node.m_points[1] - box_node.m_points[0]),
//                         glm::normalize(box_node.m_points[3] - box_node.m_points[0]),
//                         glm::normalize(box_node.m_points[4] - box_node.m_points[0])};
//    glm::vec4 extentsB(0.5f*std::sqrt(glm::gtx::norm::distance2(box_node.m_points[1], box_node.m_points[0])),
//                       0.5f*std::sqrt(glm::gtx::norm::distance2(box_node.m_points[3], box_node.m_points[0])),
//                       0.5f*std::sqrt(glm::gtx::norm::distance2(box_node.m_points[4], box_node.m_points[0])),
//                       0.0);

//    glm::vec4 translation = box_node.m_origin - box.m_origin;

//    glm::vec4 translationInA( glm::dot(translation,GToA[0]),
//                              glm::dot(translation,GToA[1]),
//                              glm::dot(translation,GToA[2]),
//                              0.0f );

//    glm::mat3x3 R;  // B to A

//    for(int i = 0; i < 3; i++) {
//        for(int k = 0; k < 3; k++)
//        {
//            R[k][i] = glm::dot(GToA[i], GToB[k]);
//        }
//    }
//    float ra, rb, t;
//    long i, k;
//    // Separatting Axis
//    // A basis vector
//    for( int i = 0; i < 3; i++ ) {
//        ra = extentsA[i];

//        rb = extentsB[0]*std::abs(R[0][i]) + extentsB[1]*std::abs(R[1][i]) + extentsB[2]*std::abs(R[2][i]);

//        t = std::abs( translationInA[i] );
//        if ( t > ra + rb ) {
//            return Outside;
//        }
//    }
//    // B basis vector
//    for( int k = 0; k < 3; k++ ) {
//        ra = extentsA[0]*std::abs(R[k][0]) + extentsA[1]*std::abs(R[k][1]) + extentsA[2]*std::abs(R[k][2]);
//        rb = extentsB[k];

//        t = std::abs( translationInA[0]*R[k][0] + translationInA[1]*R[k][1] + translationInA[2]*R[k][2] );

//        if( t > ra + rb )
//            return Outside;
//    }

//    // 9 cross products - After testing agains each obb axis

//    // L = A0 x B0
//    ra =
//    extentsA[1]*std::abs(R[0][2]) + extentsA[2]*std::abs(R[0][1]);

//    rb =
//    extentsB[1]*std::abs(R[2][0]) + extentsB[2]*std::abs(R[1][0]);

//    t =
//    std::abs( translationInA[2]*R[0][1] -
//         translationInA[1]*R[0][2] );

//    if( t > ra + rb )
//        return Outside;

//    //L = A0 x B1
//    ra =
//    extentsA[1]*std::abs(R[1][2]) + extentsA[2]*std::abs(R[1][1]);

//    rb =
//    extentsB[0]*std::abs(R[2][0]) + extentsB[2]*std::abs(R[0][0]);

//    t =
//    std::abs( translationInA[2]*R[1][1] -
//         translationInA[1]*R[1][2] );

//    if( t > ra + rb )
//        return Outside;

//    //L = A0 x B2
//    ra =
//    extentsA[1]*std::abs(R[2][2]) + extentsA[2]*std::abs(R[2][1]);

//    rb =
//    extentsB[0]*std::abs(R[1][0]) + extentsB[1]*std::abs(R[0][0]);

//    t =
//    std::abs( translationInA[2]*R[2][1] -
//         translationInA[1]*R[2][2] );

//    if( t > ra + rb )
//        return Outside;

//    //L = A1 x B0
//    ra =
//    extentsA[0]*std::abs(R[0][2]) + extentsA[2]*std::abs(R[0][0]);

//    rb =
//    extentsB[1]*std::abs(R[2][1]) + extentsB[2]*std::abs(R[1][1]);

//    t =
//    std::abs( translationInA[0]*R[0][2] -
//         translationInA[2]*R[0][0] );

//    if( t > ra + rb )
//        return Outside;

//    //L = A1 x B1
//    ra =
//    extentsA[0]*std::abs(R[1][2]) + extentsA[2]*std::abs(R[1][0]);

//    rb =
//    extentsB[0]*std::abs(R[2][1]) + extentsB[2]*std::abs(R[0][1]);

//    t =
//    std::abs( translationInA[0]*R[1][2] -
//         translationInA[2]*R[1][0] );

//    if( t > ra + rb )
//        return Outside;

//    //L = A1 x B2
//    ra =
//    extentsA[0]*std::abs(R[2][2]) + extentsA[2]*std::abs(R[2][0]);

//    rb =
//    extentsB[0]*std::abs(R[1][1]) + extentsB[1]*std::abs(R[0][1]);

//    t =
//    std::abs( translationInA[0]*R[2][2] -
//         translationInA[2]*R[2][0] );

//    if( t > ra + rb )
//        return Outside;

//    //L = A2 x B0
//    ra =
//    extentsA[0]*std::abs(R[0][1]) + extentsA[1]*std::abs(R[0][0]);

//    rb =
//    extentsB[1]*std::abs(R[2][2]) + extentsB[2]*std::abs(R[1][2]);

//    t =
//    std::abs( translationInA[1]*R[0][0] -
//         translationInA[0]*R[0][1] );

//    if( t > ra + rb )
//        return Outside;

//    //L = A2 x B1
//    ra =
//    extentsA[0]*std::abs(R[1][1]) + extentsA[1]*std::abs(R[1][0]);

//    rb =
//    extentsB[0] *std::abs(R[2][2]) + extentsB[2]*std::abs(R[0][2]);

//    t =
//    std::abs( translationInA[1]*R[1][0] -
//         translationInA[0]*R[1][1] );

//    if( t > ra + rb )
//        return Outside;

//    //L = A2 x B2
//    ra =
//    extentsA[0]*std::abs(R[2][1]) + extentsA[1]*std::abs(R[2][0]);

//    rb =
//    extentsB[0]*std::abs(R[1][2]) + extentsB[1]*std::abs(R[0][2]);

//    t =
//    std::abs( translationInA[1]*R[2][0] -
//         translationInA[0]*R[2][1] );

//    if( t > ra + rb )
//        return Outside;

//    /*no separating axis found,
//     the two boxes overlap */

//    return Partially;
}

CollisionState Collision::cube_cube_collision(Box &cube1, Box &cube2)
{
    int inside_count = 0;
    glm::vec4 local_p;
    Box *temp_box;
    temp_box = const_cast<Box*>(&cube1);
    glm::mat4 world_to_local = temp_box->m_frame->world_to_local();
    for (int i = 0; i < 8; i++) {
        local_p =  world_to_local * cube2.m_points[i];
        if(cube1.inside(local_p))
            inside_count++;
    }
    if (inside_count == 8)
        return Inside;
    if (inside_count > 0)
        return Partially;

    inside_count = 0;
    for (int i = 0; i < 8; i++){
        if(cube2.inside(cube1.m_points[i]))
            inside_count++;
    }
    if (inside_count > 0)
        return Partially;
    return Outside;
}

/*
// A simple sphere - cube collision
CollisionState Collision::sphere_cube_collision(const Sphere &sphere, const Box &cube)
{
    // TODO: use the object's frame to calculate based on the object's scale, orientation etc.
    // TODO: Change all the attributes with type Point3f to glm::
    glm::vec4 distance = sphere.m_origin - cube.m_origin;

    // Local axis relative to the box
    // X Local axis
    //Point3f xLA_temp = cube.m_points[1] - cube.m_points[0];
    glm::vec4 xLA = cube.m_points[1] - cube.m_points[0];
    // Maximum collision distance is half of the size by this local axis + the sphere radius
    float maxXLLength = std::sqrt(glm::gtx::norm::length2(xLA))*0.5 + sphere.m_radius;
    // Same for the other axis
    //Point3f yLA_temp = cube.m_points[3] - cube.m_points[0];
    glm::vec4 yLA = cube.m_points[3] - cube.m_points[0];
    float maxYLLength = std::sqrt(glm::gtx::norm::length2(yLA))*0.5 + sphere.m_radius;
    //Point3f zLA_temp = cube.m_points[4] - cube.m_points[0];
    glm::vec4 zLA = cube.m_points[4] - cube.m_points[0];
    float maxZLLength = std::sqrt(glm::gtx::norm::length2(zLA))*0.5 + sphere.m_radius;

    for( int i = 0; i <8; i++ ) {
        //qDebug() << "P[" << i <<"] : ( " << cube.m_points[i].x << cube.m_points[i].y << cube.m_points[i].z << ") ";
    }

    // For each one of the box's axis - Separating axis theorem

    //qDebug() << "Max Local X: " << maxXLLength << " - Distance in this axe: " << std::abs(glm::dot(distance, glm::normalize(xLA)));
    if( maxXLLength < std::abs(glm::dot(distance, glm::normalize(xLA))) ) {
        //qDebug() << "OUTSIDE!!!!!!!!!!!!!!!!!";
        return Outside;
    }
    //qDebug() << "Max Local Y: " << maxYLLength << " - Distance in this axe: " << std::abs(glm::dot(distance, glm::normalize(yLA)));
    if( maxYLLength < std::abs(glm::dot(distance, glm::normalize(yLA))) ) {
        //qDebug() << "OUTSIDE!!!!!!!!!!!!!!!!!";
        return Outside;
    }
    //qDebug() << "Max Local Z: " << maxZLLength << " - Distance in this axe: " << std::abs(glm::dot(distance, glm::normalize(zLA)));
    if( maxZLLength < std::abs(glm::dot(distance, glm::normalize(zLA))) ) {
        //qDebug() << "OUTSIDE!!!!!!!!!!!!!!!!!";
        return Outside;
    }

    int inside_pts = 0;
    for (int i = 0; i < 8; i++ ) {
        if (sphere.inside(cube.m_points[i])) {
            inside_pts++;
        }
    }

    //qDebug() << "Number of points inside the sphere: " << inside_pts;
    if ( inside_pts == 8 ) {
        return Inside;  // Box Inside sphere
    } else {
        return Partially;   // Sphere Partially or fully inside box
    }
}*/
